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Abstract 

In this review we first introduce the general methods to calculate the entanglement entropy for free 
fields, within the Euclidean and the real time formalisms. Then we describe the particular examples 
which have been worked out explicitly in two, three and more dimensions. 
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1 Introduction 



The entanglement entropy associated to a region V of the space in quantum field theory (QFT), is a 
quantity with manifold interests, ranging from the area of quantum computing and the study of quantum 
systems in terms of their information content [Ul^], to the physics of black holes [31I11[S] and the holographic 
principle [B] , passing through investigations with interest intrinsic to condensed matter [71 [HI E] and QFT 
itself [IU1[II1[I1]. 

The problem of an observer who has access only to a subset of the complete set of observables associated 
to a quantum system gives us the precise scenario to introduce the concept of the entanglement entropy. 
The trace over the degrees of freedom localized on a region which is not accessible to the observer, results 
in a reduced density matrix. In particular, if the state of the quantum system is the vacuum, the local 
density matrix pv reduced to a region V of the space is 

pv - tr_y |0) (0| , (1) 

where the trace is over the degrees of freedom lying outside V. Then, the entanglement entropy is the von 
Neumann entropy associated to pv 

S{V)=-tv{pv\ogpv). (2) 

The free quantum field theories play here, as is often the case, the role of the simplest models available 
to the research. However, even in the free case exact results on the entanglement entropy are scarce, and 
require solving difficult questions about the spectrum of certain integral operators, or equivalent questions 
expressed as functional integrals. Powerful numerical methods on the lattice are known which reduce the 
calculation of S(V) to a diagonalization of few ^ mx m matrices, where m is the number of lattice points 
in V. This is in contrast to the ^ q™' x q"^ matrix which is at least expected in a generic model with q 
states per site. This important reduction in computing time allows one to take a close look at the behavior 
of the entropy function in the free case, at least for low dimension space-times. 

The difference between free and interacting QFT models is not evident in terms of the entanglement 
entropy alone, and has not been elucidated yet in these terms. This makes the study of the free models 
more interesting because it may reveal features which are common to all QFT. In particular, the ultraviolet 
behavior is one of these features. In d spatial dimensions we should have for any QFT 

S{V) = gd-i[dV] e-('^-i) + ...+gi[dV] e'^ + go[dV] log(e) + So{V) , (3) 

where 5*0(1^) is a finite part, e is a short distance cutoff, and the gi are local and extensive functions on 
the boundary dV, which are homogeneous of degree i. The leading divergent term coefficient gd-i[dV] is 
then proportional to the {d — 1) power of the size of V, and this is usually referred to as the area law for 
the entanglement entropy. 

This area law is a consequence of the large number of degrees of freedom at high energies present in the 
QFT which induce entanglement across the boundary dV. The divergent terms are essentially produced 
by high energy vacuum fluctuations and thus they are the same for all finite energy density states, such 
as a thermal stateQ 

These terms proportional to gi for i > are not physical within QFT since they are not related to 
quantities well defined in the continuum. On the contrary, the coefficient go of the log term is expected to 
be universal (in this review we use the word universal in the sense of independence of the regularization 
prescription or of the microscopic model used to obtain the continuum QFT at distances large with respect 
to the cutoff). There are also many universal terms which are included in 5*0, but have to be extracted 
taking into account that the definition of 5*0 is affected by finite changes in the cutoff. 

Among the universal quantities related to the entropy an important role is played by the mutual 
information I{A,B) between two non intersecting regions A and B T], (see figure [T]) 

I{A,B)^SiA) + S{B)-S{AuB). (4) 

^ There is a different meaning of area law wliicli is used is the literature of discrete systems which in fact concerns the 
behavior of the entropy in the large volume limit (rather than a fixed volume in the continuum limit as in (|3j) for specific 
states such as the ground state, or a state with chemical potential (for a review and a complete list of references see 1131 V 
It is known that for fermions in a lattice and finite chemical potential the entanglement entropy grows as the area with a 
multiplicative logarithmic correction 1141 . However, this does not affect the formula (|3]l since these terms induced by the 
chemical potential are finite and independent of the cutoff in the continuum limit, that is, they are included in the finite part 
So- The same holds for the finite volume increasing terms in a thermal state, which should also satisfy an area law in the 
sense of ((Sj. 
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Figure 1 : The figure on the left shows the term 5*0 in the entropy of circles in two dimensions as a function 
of the radius R, with a square lattice regularization. 5*0 is obtained by fitting the data as S{R) = ciR + Sq, 
and subtracting the term ciR. It is apparent that 5*0 does not converge to any definite value. The figure 
on the right shows the mutual information I(A, B) between two circles A and B of radius R on the lattice, 
separated by a distance to each other which we have chosen to be also R. The calculation shows the 
convergence of the mutual information in the continuum limit i? — > cx) (with an error of ~ 2% already for 
R ^ 30). In both figures we have considered a massless scalar field and R is measured in lattice units. 



The universal character is due to the cancellation of the boundary terms. Another useful dimensionless 
universal quantity defined in two dimensions and for an interval of size L is given by 

c(L).L«:l. (5) 

The entropic c-function c{L) is always positive and decreasing, and plays the role of the Zamolodchikov's 
c-function [15j in the entanglement entropy c-theorem |12| . It contains all the universal information present 
on the entropy for an interval, since S{L) follows from c(L) by integration, except for an arbitrary constant. 
Accordingly, in this review we will present the results for the entanglement entropy for single intervals in 
two dimensions in terms of c(L). 

Still related to the local density matrix, there is a family of functions called the alpha or Renyi entropies, 

^n(V^) = -^logtrp^, (6) 

1 — n 

which share some of the entanglement entropy properties. These are often easier to compute, and one has 

lim5„(y) = (7) 

For an interval in two dimensions a quantity analogous to the entropic c-function ([5]) can be defined as 

c„(L) = L^%^, limc„(L) = c(i). (8) 

This review is focused on vacuum entanglement entropy in relativistic free QFT and its associated 
universal quantities. In what follows we describe two general methods to calculate the geometric entropy 
for free fields. We then give examples where these methods are applied in two, three and more dimensions 
for scalar and Dirac fields. Free gauge fields have not been sufficiently studied in the literature yet. They 
have an additional technical problem since the divergent terms in ([3]) are non gauge invariant. 

Most of the material is already present in the literature, but we have also included some new results, 
and reformulated earlier work. There are several topics related to entanglement for free fields discussed 
in the literature but not belonging to our main subject which will not be covered here. These topics 
include entanglement entropy in presence of temperature, chemical potential, or excited states [TOl [Ml 
[ini [TT], manifolds with boundary [T51 [121 [201 [H] , different (Newton- Wigner) localization schemes |22| . 
entanglement in curved space-time [231 124j , and in non-Lorentz covariant QFT (quantum Lifshitz fixed 
points) [53]. There is also a large body of work done in the context of lattice models, which extends out 
of the scope of this review. There are several excellent review papers in this area [51 [T51 [^ . 
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2 The entanglement entropy in free QFT 



2.1 Euclidean time method 

In the Euclidean approach one starts from the representation of the vacuum state in terms of a path 
integral. For definiteness consider a scalar field (/)(af, t), and take the basis formed by eigenvectors of this 
field operator at time t = as (j){x,0) \a) = a{x) \a), where a is any real function on the space. The 
vacuum wave functional writes [27j 

= (0|a) = / e-Si,(0) (g) 

<p{x, — oo)—0 

In order to select the vacuum state, the functional integral is over the lower half space and with Euclidean 
time. SEi(t>) is the Euclidean action and N^^/"^ is a normalization factor. The vacuum density matrix in 
this basis is p{a,a') — (a|0) (Oja') = <i>(Q!)*<i>(Q!'). In order to trace over degrees of freedom in —V, the set 
complementary to V, one considers functions a = /3®Q;y, a' = [3®a'y, which coincide (are equal to [3) on 
— V, and sum over all possible functions (3. Using the representation ^ this construction of the reduced 
density matrix amounts to take two copies of the half space, glue them on —V (see figure^]), and take the 
functional integral in this space [55], 

/l'(t>{x,0'^)=av(x),x£V 
DP ^{(3®avTHP®oi'y)=N-^ D(f)e-^^^'^K (10) 

J cj>(x,0-)=a'y(x),xeV 

The arguments of the density matrix are the boundary conditions of the path integral on both sides of the 
cut on V . 

For fermions a similar construction holds, but there is an important difference. Since the fields at 
equal time anticommute, the functional integral and the boundary conditions are in terms of Grassmann 
variables. In order to represent the trace in the functional integral we have to sum over the antidiagonal 
elements [29], that is 

pv{av,a'v) = j D(3 'i>{-l3®avTHP®a'y), (11) 

where now all variables involved are Grassmann valued. This extra minus sign is well known in the 
calculation of the thermal partition function for fermions, which requires antiperiodic boundary conditions. 
Here, this implies that when we glue the two copies of the half space along —V we have to take boundary 
conditions related by a minus sign. This cut at — can be eliminated by changing integration variables 
^ —ip in the upper half plane. The net result is that there is an additional minus sign in the boundary 
condition for the density matrix in the upper side of the cut over V in figure [U 

l'ip{x.O''')=-av{S),x<£V 

pv{av,a'y) = Di^ e'^''^'^ . (12) 

J ■4i{x,0-)=a'y{x),x£V 

The representation of the traces trpy involved in ([S]) by a functional integral is then realized by the 
replication method which consists on taking n copies of the Euclidean plane cut along V ^ and sewing 
together the upper side of the cut in the fc*'* copy with the lower one of the (fc + 1)*'* copy, for fc = 1, n 
[H [To] . Here the copy n + 1 coincides with the first one (see figure [3]) . The resulting space is a n-sheeted 
d + 1 dimensional Euclidean space with conical singularities of angle 27rn located at the boundary dV of 
the set V . For fermions, the final trace (after the power has been taken) also involves a sign change 
in the boundary conditions. We can change variables in order to have all gluing with no additional signs 
except for the last one, where the field on the upper side of the n*'' cut is equated to (—1)"+^ times the 
field on the lower side of the first cut [30|. Here, one minus sign counts for each copy of the density matrix 
plus one extra minus sign for the traco 



^For the special case of the Rindler space, corresponding to V half a spatial hyperplane, the density matrix is proportional 
to e~^'"^ , with K the boost operator which keeps the Rindler wedge fixed }31J . In the Euclidean version, the boost operator 
corresponds to a rotation operator and p" = g— ^-n-nif ^ rotation operator of angle 2-Kn. When translated into the path 
integral formalism this provides a different look to the origin of the signs on the fermion boundary conditions in this case 
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Figure 2: The path integral on the lower half Euclidean space gives the vacuum wave functional. The re- 
duced density matrix pv is obtained by gluing two copies of this half space along the set — V complementary 
to V. 



Finally we have 



Sn(V) 



Z{n) 
Z{1)" ' 

logZ(n) - nlogZ(l) 
1 — n 



(13) 
(14) 



where Z(n) is the functional integral on the n-sheeted manifold, and we have used the normalization factor 
N = Z{1) in order to have tipv = 1- Eq. I|14[) gives a representation of the Renyi entropies for integer n. 
The entanglement entropy follows by analytic continuation of Sn down to n — 1 

2.1.1 Diagonalization in replica space 

In general, calculating the integrals Z{n) explicitly is a very difficult problem since we have to deal with 
a non trivial manifold resulting from the replication method. Fortunately, in the case of free fields, a 
simplification follows by mapping the n-sheeted problem to an equivalent one in which one deals with n 
decoupled and multivalued free fields [30 . For that, one introduces a vector field $ living on a single d+1 
dimensional space, whose components are the values of the fields in the different copies. 




(15) 



where 0/(a;) is the field on the I copy. Note that in this way the space is simply connected but the 
singularities at the boundaries of V are still there since the vector $ is not singled valued. In fact, crossing 
V from above or from below, the field gets multiplied by a matrix T or respectively, where 



/ 



V (±1 



in+1 



1 

/ 



(16) 



The upper sign in this equation corresponds to the bosonic case and the lower one to fermions. 

This matrix has eigenvalues e*"^'^, with k — 0,..., {n — 1) (the n^^ roots of 1) in the scalar case, and 
with k ^ -{n- l)/2,-(n- l)/2 + 1,..., (n - l)/2 (the n'^ roots of (-1)"+^), in the fermionic 
one. Then, changing basis by a unitary transformation in the replica space, we can diagonalize T, and the 
problem is reduced to n decoupled fields 0^ living on a single d + I dimensional space. At this point is 
essential the free character of the action, and we also need to deal with complex fields. This last requirement 
is not a limitation since in order to compute the entanglement entropy we can double the number of real 
fields and then divide the final result by two. The fields which diagonalize T are multivalued and defined 
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Figure 3: tvpy is given by the path integral on a n-sheeted space formed by sewing the rephcated Euchdean 
spaces with a cut along V. This is equivalent for free fields to n decoupled multivalued fields living in a 
single space. These fields get multiplies by particular phase factors when crossing the cut. 



on the Euclidean d + 1 dimensional space with boundary conditions imposed on the d dimensional set V 
given by 

<^fe(x,0+)-e^^<^fe(f,0-) , xeV. (17) 

Here (pkix, 0*) are the limits of the field as the variable approaches V from each of its two opposite sides 
in d + 1 dimensions. 

In this formulation we have for scalars 

n— 1 

Sn{V) = E log^[e''"'^"] , (18) 

and for fermions 



S„iV) = -^ V logZ[e'2.fc/„j, (19) 
1 — n ^ — ' 

fc=-(n-l)/2 

where Z[e'^'^''] is the partition function corresponding to a field which acquires a phase e'^'^° when the 
variable crosses V (figure [3]), divided by Z{1). 

To further specify these partition functions one has to take into account that the fields must have a 
specific asymptotic behavior as the variable approaches the singularity in order to have finite action. This 
requirement will be different according the action is quadratic (scalars) or linear (fermions) in derivatives. 
For smooth V boundaries we have 

(f>{x) - r'^ , 7 > (scalars) , (20) 

ip{x) ~ r'^, 7>-l/2 (fermions), (21) 

in the limit of short r, where r is the distance between x and dV . 

Note that (fTSl) and (|19p give the same formula for the limit n — > cx), 

So.iV) = -^ J d9 logZ[e'Y (22) 



2.1.2 Analytic continuation 

The evaluation of the entropy is still limited by the difficulty in doing the analytic continuation of S„ for 
non-integer n, and then the limit ([7]). The analytic continuation of a function defined on the integers is 
not unique unless some further information is provided. One such requirement would be the condition 
Sn < ce'^l"! for Re(n) > 1/2, which is the hypothesis of the Carlson's theorem [33l[34]. This holds for 
finite dimensional density matrices, but for universal terms in 5„ in a continuum theory no general result 
is known. For a discussion around this point see [331 135j . 

The analytic continuation for the free fields which satisfies Carlson's criterion can be obtained in a 
natural way writing the sums in Sn as a contour integrals [33]. We will check in section 2.3 that the result 
coincides with the entanglement entropy obtained by the real time approach. 
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Figure 4: The integration path used in eqs. and (|27p to define the analytic continuation of trp" for 
non-integer n. 

We treat first the scalar case. Eq. (fTS]) can be written 

' ^^^du^^^^,^-^Jdu--^Xo,Z[u], (23) 



" 2m(l-n)^J u-e*2'r'=/" 2m(l - n) J u" - 1 

where the contour of integration encircles the roots of the unit. The function logZ[M] is defined as the 
unique analytic continuation of the function log Z[e*^'^''] on the unit circle. Since we need to avoid the 
negative real axis for the formula (|23p to be analytic for non integer n, we choose the integration contour 
as in figure[4l Thus, we are assuming logZ[u] is analytic on the negative real axes and |logZ[M]| / \u\^ 
as u — oo for some x < 1 (we will show this holds for the universal terms in log.Z[u] in section 2.3). In 
the limit n — > 1 we have 

- i /"« (;7^T^ + (Sf + - ■ <^'') 

At M = 1, logZ[ii] vanishes. We also have the identity Z[e'^'^"] = Z[e~*^'^°], due to the EucHdean time 
reflection symmetry. This leads to 

Z[u] = Z[l/u]. (25) 

Then, logZ[u] should be at least of order {u — 1)^ near u = 1. In consequence, the first term within the 
brackets does not have singularities inside the contour and does not contribute. Thus, the only singularity 
is the cut of the logarithm and we obtain 

where we have used the property (|25p . 

The fermion case can be treated similarly. The sum of residues in this case is over the roots of (— 1)""*"^, 
giving 

For the last integral in ?j = —u we choose the contour in figured (the same as for scalars). Again, using 
property the limit n ^ 1 gives 

5(^,.2fdA;^, (28, 

where we assume |logZ[A] | /A"^ ^ as A ^ oo for some x < \. In concrete examples, the result for 
log Z[e*^'^"] automatically gives a formula for the analytic extension log Z[u]. Then the entropy is calculated 
by (Ell) or (EHl)- 
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2.1.3 Methods for calculating the partition function 

One of the most powerful methods for computing the partition function of quadratic actions is the heat 
kernel method. To introduce it, consider a free scalar, and note that 

W = ~log(Z) = ilogdct(m2 - V^) = itrlog(m2 - V^). (29) 

The partition function then reduces to an spectral function of a differential operator. The heat kernel is 
defined as K{x,y,t) = {x\ e*^ \y), and its trace writes 

C{t) ^tre*^^ ^ J dxK{x,x,t) . (30) 

The free energy W can be written in terms of this spectral function as 

W^-^J^ dte-^-''-C{t), (31) 

where e is a cutoff. The advantage of this equation if that, in contrast with (|29p . it is written as a function 
of the trace of an operator which satisfies a local heat equation 

dK 

— = V^K, K{x,y,0)^S{x~y). (32) 

This means that systematic expansions for short time t of the trace of the heat kernel can be developed, 
leading to an expansion of the form (for a review see [37]) 

C(t) = ^i('=-^)/V. (33) 

fc>0 

Here D = d + 1 is the manifold dimension, and d is the dimension of space. In general, the heat kernel 
coefficients Ofc are integrals of local quantities depending on the different background tensors. This may 
include a smooth metric, interactions with external sources or fields, or a smooth boundary and boundary 
conditions. An expansion is known for different spins up to the first few k in the general case. 

However, the application of the heat kernel expansions to the entanglement entropy in flat space vacuum 
is made difficult because of two reasons. The first one is that the short t expansion is an ultraviolet one. 
Thus, the non divergent contributions, which represent genuine physical terms in the entanglement entropy, 
are the ones proportional to the coefficients Uk with k > d + 1. The case k = d + 1 gives a logarithmically 
divergent contribution to the entanglement entropy, with universal coefficient which is proportional to 
the conformal anomaly (see section 3.3.3). Obtaining these large-index coefficients becomes increasingly 
difficult for higher dimensions. 

The second reason is that the manifold has conical singularities along the boundary of the region V, 
and the standard expansions cannot be applied directly to this case. Expansions in presence of conical 
singularities in the limit of small deficit angle, with bulk curvature, and smooth dV having zero extrinsic 
curvature, have been developed in [55] (see also [35]). They where applied to obtain logarithmic corrections 
of the black hole entropy in different scenarios [531 [35] . Recently, the contribution of the extrinsic curvature 
was found in four dimensions, and used to calculate the logarithmic term of the vacuum entanglement 
entropy in Minkowski space with smooth dV 4^. We will review this result in section 3.3.3. 

Contributions from a non smooth boundary become intractable with the heat kernel method. We will 
encounter this type of term in section 3.2. 

Another possible strategy to find the partition function consists in studying the associated Green 
function G = (— + w?)~^ on the manifold relevant to the problem. These are related by the identity 

^logZ = -itrG. (34) 

However, there are no general methods to calculate G for this type of manifolds with a codimension one 
cut on a finite region, and a case by case analysis seems unavoidable. In some of the examples we review 
in section 3, we have used a method originally introduced by Myers in [41] to deal with Neumann and 
Dirichlet boundary conditions on a finite cut in two dimensions. It essentially consists in exploiting the 
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symmetries of the Helmholtz equation, even in the presence of the boundary conditions, by analyzing the 
singular behavior of G(r, r') at the end points of the cut, which effectively control the solutions of the 
source free equation. With this tool we could calculate the entropy for a single interval for massive free 
fields (sections 3.1.1 and 3.1.2) and the coefficient of the term logarithmic in the cutoff for a plane angular 
sector in two dimensions (section 3.2.1 and 3.2.2). 

In d = 1 other techniques are also available. In particular the bosonization technique can be used to 
evaluate Z{n) for the Dirac field (section 3.1.5). One expresses the fermionic current in terms of a dual 
scalar field cj) as -^e^'^d^4). The resulting dual theory is the Sine Gordon theory, and the functional 

Z{n) is given in terms of a sum of correlators of local exponential operators [3D]. This is a particular case 
of the general one in one dimensional space, since Z{'n) can always be expressed as a correlator of twisting 
operators [33] (see also [H]). These are non local relative to the ordinary fields, and effectively impose the 
boundary conditions. The correlators of twisting operators for general integrable massive quantum field 
theories have been studied from a S matrix approach. In particular cases, an expansion of the correlators 
in form factors is known, and they can be evaluated as a sum over intermediate multi-particle states. These 
expansions have been exploited in relation with entanglement entropy in [191 1331 135j . 



2.2 Real time approach 

In the real time approach one aims to compute directly the reduced density matrix corresponding to the 
global vacuum state in terms of correlators. One starts with a discrete version of the quantum field theory 
and eventually takes the continuum limit. The first calculations of the entanglement entropy where made 
in this way [3] , several years before the Euclidean approach was developed. This method has been mainly 
applied to numerical calculations in the lattice. Even if it has not been exploited as much as the Euclidean 
approach for analytic calculations, we think it is better suited to several generalizations. For example, it 
can be used in calculations involving spatial sets which are not contained in a single spatial hyperplane 
in Minkowski space, for which the Euclidean approach becomes inapplicable. The interactions can be 
included in a straightforward way, at least at the perturbative level \AS . Also, the case of some states 
different from the vacuum seems to be tractable with this method. However, these applications have not 
been fully developed in the literature yet. 

First, let us describe briefly how this approach was formulated originally. For a system of harmonic 
oscillators with Hamiltonian 

N N 

H^-Y.i'l + i^Y.'^^^^o'Pj (35) 
the ground state wave function can be easily obtained by diagonalization 

M.o(0)= (^det-^)^-^^"^^ (36) 

where W — \Ik. The corresponding density matrix in coordinate representation is 



p(0, <iJ) = *S(</')*o(0') = ,/det— e-5^"^*e-5^'"^^' . (37) 



TT 



Then, one selects the variables 0'-*^ and (\)^°^ corresponding to the degrees of freedom inside and outside 
V . The trace over the 0*^°^ can be done by integration of the Gaussian in (|37)) . Writing 

M/=f5« (38) 
in the base where the inside variables g'*-' are the first ones, one arrives at 



py(x,x') = ydet^— e-^" "e-^"' -'e3("+"') M^+^O ^ (39) 

where 

X = (40) 

A = w^(-;/V(,<,)I¥(;lw^(„)I^-//^ (41) 
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The entropy follows by writing the density matrix as one for non coupled degrees of freedom by making a 
linear transformation of coordinates. Finally we have 




siv) = tr I iog( \-_^f^+;^ ) - v:-;^: i ■ (42) 

This calculation has first been done in [5]. It was later rediscovered in [5] where it was applied to free 
quantum fields discretized in the radial direction in polar coordinates. For fermions a similar calculation, 
has been done in [9]. 

In the following, we use a method which starts with a different perspective, and gives equivalent 
results in a way more suitable to generalizations. By its very definition, the reduced density matrix pv 
corresponding to the region V, is the state acting on the local algebra of operators in V which leads to the 
same expectation values as the global vacuum state, 

(Ov) = tr(pvOy) , (43) 

for any operator Oy localized inside a V. If we take the point of view of the Wightman theorem [31], where 
one defines a quantum field theory in terms of the correlators, this equation suggests that the knowledge 
of all the correlators inside V must be enough to determine the density matrix pv- This is specially 
simple in the free case, where the Wick theorem means that all the correlators are reduced to the two 
point function. This fact was used by Peschel to give an expression for pv in terms of correlators for free 
Boson and Fermion discrete systems |45j (see also j9j). Similar expressions for the local density matrices 
in terms of correlators where obtained long time ago by Araki |46j . His work is framed in a more general 
and mathematically oriented context, classifying the possible states satisfying the Wick theorem (called 
quasifree states) for algebras obeying the canonical commutation and anti-commutation relations. 

2.2.1 Bosons 

The local Hermitian variables (f>i and nj (coordinate and conjugate momentum) obey the canonical com- 
mutation relations 

[(f>i, TTj] = iSij , [(l>i,(f>j] = [tTi, TTj] = . (44) 

Take the two point correlators inside V to be 

{(l)t(f)j) = Xij , (TTiTTj) = Pij , (45) 

2 

" TT,,) = (TTj^i)* = . (46) 



■Jin J 



This last equation can be generalized to have {4>iT^j) -f {'^j4'i) 7^ 0, but (j46p is all we need for the vacuum 
state. The equations in (|45[) imply the matrices X and P are real Hermitian and positive. The positivity 
of (((/); -|- iXik'^k){4'm — 'i'Kns''^s)) for arbitrary constants Xik implies that 

X.P > i , (47) 

in matrix sense, what means that the eigenvalues of X.P are greater than 1/4. 

We assume all other non zero correlators are obtained from these by the Wick's theorem 

(T 

where the sum is over all the index permutations ti, the fi can be any of the field or momentum variables, 
and O is an ordering prescription, for example, ordering the products inside the expectation values with 
the field variables at the left and the momentum variables on the right. Once this equation holds for a 
specific ordering automatically holds for any other ordering. 

Let us consider now general independent creation and annihilation operators a;, a\, with [a^, aj] = 5y , 
which are expressed as linear combinations of the and tTj, i, j G V , 

fpi = a*ja^^ + aijUj , (49) 
TTi = —i/3*ja^j + iPijUj . (50) 
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The commutation relations between coordinates and momentum give 



a*/3^ + a/3^ = -1. 

At this point we note that an ansatz for the density matrix of the form [351 [37] 



pv = Ke 



-H 



Ke 



Seiajai 



(51) 



(52) 



where the normahzation constant K — n;(l — e~''), leads automatically to the Wick property for the corre- 
lators {Ofi^fi^...fi.^^) — iv{pvOfi^fi^...fi^^^) HE]- Here Ti. is called the modular Hamiltonian corresponding 
to the region V |49j, and ([5^ proposes it quadratic in terms of the creation and annihilation operators. 
The Wick property holds independently for every mode, as can be shown by direct computation, 



tr 



i!tr 



(1 - e-^-)e"''°'°'aIazl " = n!(e"' - 1)"" 



(53) 



From linearity, the eq. (|48p follows. Then, with this expression for pv we calculate tT[pv4>iT^j) = (i/2)(5y , 
tr{pv4>i4>j) = Xij, and tr(pv'7ri7rj) = Pij. This gives 



1 

2 ' 



/3*n/3^ + ^(n + l)/3^ =P, 



where n is the diagonal matrix of the expectation value of the occupation number 



nkk 



(54) 

(55) 
(56) 



(57) 



These equations give that a — aiU and (3 — (3iU , with U unitary and diagonal, and j3i and ai real. 
The matrix of phases U can be reabsorbed in the definition of the a,;, so we set U = \. Then we have 

^'»T^^l and 



0' 



a-(2n+\fa-^ = XP . 
4 



This last equation gives the spectra of the density matrix in terms of the spectrum of XP^ 

(1/2) coth(efc/2) = , 



(58) 



(59) 



where Vk are the (positive) eigenvalues of C = \J XP . This means that, in the bosonic case, the modular 
Hamiltonian is definite positive. 

Inverting the relations ((33) and ([50]) and replacing in ([5^ we can write the density matrix as 



PV 



where 



M 



N 



^ -IT -1 

-a ta 



^2C1°S 



C 



aea 



2C^°^ 



C 



X. 



(60) 

(61) 
(62) 



with e the diagonal matrix of the e^. The entropy is given by 



5 = ^(-log(l 



tr ((C + 1/2) log(C + 1/2) - (C - 1/2) log(C - 1/2)) , (63) 



1 - e-" 

which is positive thanks to C > 1/2, eq. (j47p . We also have 

log(trp«) = -tr [log ((C + 1/2)" - {C - 1/2)")] 



(64) 
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For a quadratic Hamiltonian for bosonic degrees of freedom with the form (|35p , H — ^ '"'I + 
i 4'iKij4'ji the vacuum (ground state) correlators are given by 

X,, = {^,cj,,)^^{K-i),,, (65) 

P,^ = (^,7r,) = i(if^),,. (66) 

For the global state we have X.P — 1/4, which has zero entropy, corresponding to a pure state. 

A straightforward calculation shows that in this case C = VXP = 1/2(1 - A)-i/2 (where A IS given 
by (|4T|l ). and the expression ((63l) for the entropy is equivalent to (|42|l obtained with the Hamiltonian 
approach. However, this formulation stresses the fact that we need only the correlators inside V, which 
is not apparent in the eq. (|4ip . involving matrix elements in the exterior region. In some situations, 
the field correlators are known, and this knowledge may lead to a great reduction in computing time in 
numerical simulations, avoiding the use of an infrared cutoff. Besides, it makes more transparent the range 
of applicability of the formalism, which extends to all states satisfying the Wick theorem (Gaussian states). 
This includes states in thermal equilibrium and with chemical potential, or some vacuum states in curved 
space-time. 

In the continuum limit C converges to an integral kernel. However, in this limit is not obvious how this 
depends on the correlators {(t){0)(f){x)) and (7r(0)7r(a;)), since these do not define good operator kernels. 



2.2.2 Fermions 

The local creation and annihilation operators "ipj, ipj satisfy the anticommutation relations {ipi, V^j} = Sij. 
Let the two point correlators be 

^.v-l) = aj, ^ - c,, , (67) 

(V'.Vj) = (Vl^])=0. (68) 

We also assume the Wick theorem holds and all non-zero multipoint correlators are obtained from the two 
point functions in the standard way 

k 

(v'n...^..^L..V'I) - (-1)"^ n {^^.^L^) ' (69) 

cr q—1 

where the sum is over all permutations of the set of indices ji, ■■■,jk, and is the permutation signature. 
The case when the correlators of Eq. ([55]) are different from zero can also be treated, but we will not 
need this generality here. The hermitian matrix Cij of correlators inside the region V must determine 
completely the local state pv- From ([57)1 Cij and 5ij — Cji are positive, and thus have eigenvalues in the 
interval [0, 1]. For V the total space, and when the global state is pure, C is a projector, having eigenvalues 
or 1. 

Analogously to the bosonic case, the correlators computed with the help of the reduced density matrix 
Pv and eq. (|43)) satisfy the Wick property if we have [48] 



PV 



= Ke-^ = Ke- ■'^^^ (70) 



Since pv must be Hermitian H, is also Hermitian. We can then diagonalize the exponent by the Bogoliuvov 
transformation di — Uim'4'm, with U unitary in order to maintain the anticommutation relation, {di, dj} = 
5ij. One then chooses U such that U.H.W — {e;} is a diagonal matrix where the are the eigenvalues of 
H. One has 

and from here we specify the constant K = det(l -f e^^)^^. 
The relation between H and C follows from 



K tr 



(^e-^v^.^^!^^^^^^^^=C.^, (72) 
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Using the diagonalization ((7T|l one obtains that the eigenvalues e/ of H and the eigenvalues vi of C are 
related by 

e^'=^, (73) 

or, in matrix notation, 

H = -log(C-i-l) . (74) 

With vi e (0, 1) the one particle Hamiltonian Ti. is not defined positive. However, the large negative energies 
do not contribute much to the entropy, since the corresponding degrees of freedom become saturated at 
occupation number 1 (the equivalent to the Dirac sea). The large positive energies saturate at occupation 
number 0. 

The entropy and log trp" can be evaluated as a sum over independent modes 

5(F) = (^log(l + e-^') + Y^5^) =-tr((l-C) log(l - C) + C log(C)) , (75) 

logtrp'^ = S,(log(l + e""^')-'^log(l + e-^0) = tr (log((l - C)" + C")) . (76) 
The typical case where the Wick theorem holds is for a quadratic global Hamiltonian of the form 

7^ = 5]A/yVjV',, (77) 

i,3 



choosing any Fermi level and temperature, since in this case the global density matrix is an exponential 
of a quadratic form in the creation and annihilation operators. 

In quantum field theory applications one takes the vacuum (half filled) state of a Hamiltonian with 
symmetric spectrum around the origin. The explicit form of the correlator in this case is a projector 

C = 0{-M), (78) 

where 0{x) is the step function. For the whole space we then have zero entropy. The reduction to a region 
makes the eigenvalues of C lie between zero and one. In this case we also have a total symmetry between 
ipi and tpj, (a local charge conjugation symmetry) which leads to the same spectrum for C and 1 — C. 
Therefore the spectrum of energies of the modular Hamiltonian is also symmetric around 0. 
Density matrix for a Dirac field 

As an application of these formulas to quantum field theory, let us consider a free Dirac field in d + 1 
dimensions (an early treatment of the continuum limit for fcrmions can be found in [50j). The field 
satisfies the Dirac equation {id^'y^ — m) 'J = 0, with 7^ the Dirac matrices, and the canonical (equal-time) 
anticommutation relations 

[^,{x,t),-^]{v,t')'^\^^^^^5''{x~v)5,,. (79) 
We can discretize the field algebra (at i = 0) by using the smoothed field operators 

V'™ = / dx'*«'.(f,OX(f), (80) 



where an{x) is an orthonormal base of spinor functions on V ^ Jy dx'^ al^{x)am{x) = 5n,m, and the discrete 
field operators satisfy {i/'m, ■i/'Jj} — Sm^n- The results of the previous section can be directly applied to 
these set of discrete operators. After that one can recover the field language with the help of ([5U)) . 

Let then the set V be an arbitrary region contained in a surface of constant time. We have the 
expressions analogous to ([70)1 and ([71)1 

P = Ke-^, (81) 

n = f dx" dy^-^l{x,0)H,,ix,y)'<SJjiy,0) , (82) 
Jv 

H = -log(C-i-l), (83) 

where the field correlator is 

C(.^,.1= (0|v,.(f,0)*t(,^,0)|0)L,_,, = / (g^|^J=i|/e-*-). (84) 

The expression for the entropy is again ([75]) . but now C is an operator with kernel rather than 
a matrix, and (j75p requires a regularization. The more general case of an arbitrary spatial region not 
contained in a single spatial hyperplane is treated in [51j . 
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2.3 Direct relation between the Euclidean and the real time approaches 

The expression for the entropy in the Euchdean approach can be directly transformed into the real time ap- 
proach formulas ([55)1 and (|75p by evaluation of the functional integral Z[X] in terms of boundary operators. 
A similar treatment has been applied to the Casimir effect in [52j . 

Let us start with the Dirac field. The boundary condition in Z[X\ , with A ~ e«27ra^ jg taken into account 
with the modified action S[^, ^] = So[^, + SSv[^, where 

S-of*,*] = / dx-^+i *(7^9^, + m)* , (85) 



SSv - 



dx'^+i (1 - A)(5(r)xi/(f)*(x)7°*(a;) = (1 - A) / dx'^ *7°^' , (86) 

Jv 



where r is the Euclidean time coordinate, xvix) equal to one for x G V and outside, and the last 
integral over the spatial d dimensional set V. This is because this new action leads to the same classical 
solutions and boundary conditions as the original problem. It is also possible to understand the term 
5Sv as produced by an external gauge field vanishing outside V, which is pure gauge everywhere except 
at dV, and which has the effect of imposing the correct boundary conditions (see section 3.1.5). Strictly 
speaking, in a discretization of the path integral, the term in (|86p has to be understood as proportional to 
the product ^~^{x)^°'i/~{x) of the fields located on different sides of the cut. 

Writing the partition function in terms of an auxiliary Grassmann field living in V in order to linearize 
SSv we have 

Integrating first over ^ we have 

Z[X] f V^V^ e-^vd-' Ivdv'ii-H^i-}^iv})Bay)+T^Udx^i-f"i ^^Qtifx{l + iX-l)C)), (88) 



where, as in section 2.2.2, C is the Minkowskian correlator C{x,y) = {'^ (x)'^^ (y)) — {'i> {x)^ {y)^^) e inside 
V, and fx is an unimportant normalization factor which cannot change the universal terms in logZ[A]. 
The analytic extension of Z[X] for |A| 7^ 1 is done directly with formula ([55)1 . The symmetry Z[X] — Z[l/X] 
(eq. (HH)) (for the universal terms) is mapped here to the identity of the spectra of Cv € (0, 1) and 1 — Cv- 
We can choose fx = A~^/^ in order to have Z[l] = 1, and Z[X] = Z[l/X] exactly. Using (ISS)) in the formula 
(j28p obtained from the analytic continuation of Sn in section 2.1.3, and taking into account the spectral 
properties of C, we obtain the correct formula for the entropy (|75p . 

The boundary condition for the scalar partition function can be similarly implemented by adding a 
term to the first order action, 

S[(t>, Tr]= J dx'^+^ +71*^^ {nn* +V(t>V(t)* +m^(t)(j)*)^+ J dx''((l - A)(/.7r* + (1 - A-i)^* tt) . (89) 

Here A is again a phase factor. Writing the corresponding partition function in terms of auxiliary fields 
living on we get 

We can now integrate over tt and (j) first, giving 

z[x] = f vc^v^^vi^vi^ exp ( f dxUy''{CA^){^{x)^*{y))EUy) + C2{x){HxW{y))EUy) (9i) 



+ {C{x){7T{x)cf,*{y))EUy)+i^-C-)+ I dx" (^ + -^4^)) = dct-^ (fx (\ ~ f,^ , ^}l c^ 



1-A l-A-17 V V (1 + ^) 

where C = ^/XP, with X and P the operators with kernel given respectively by the Minkowskian corre- 
lators {(j){x)(jf {y)) = {(j){x)(j)* {y)) E and {tt{x)tt* {y)) = -{t:{x)'k* {y)) e restricted to V (see section 2.2.1). 
We have also used {Tr{x)(j>* {y))E = {<t>{x)'n:* {y)) e = -l/2S{x - y). When ([HT]) is inserted in it gives 
the correct formula for the entropy (j63p , obtained with the real time formalism. 
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It follows from (l88l) and (j9ip and the spectral properties of the operator C in each case, that the 
universal terms in logZ[A] (disregarding a global constant independent of C) satisfy the conditions for 
|A| — *■ oo which were assumed in section 2.1.3. 

We note that ([57)1 and ([50)1 give definitions for the partition functions log Z[\] which extend analytically 
out of the range |A| ~ 1. The induced boundary conditions in this case are given by a factor A for the 
fields when crossing the cut, but a factor A"-'^ for the conjugate momentum. 



3 Exact results 

3.1 One spatial dimension 

3.1.1 Single interval for a massive scalar field 

In this section we review the calculation of the universal part in S'„ and S, for a massive scalar field in a 
single interval [L2,Li] presented in |53| . We use the Euclidean time method and calculate the partition 
function Z[X] for a complex scalar from the Green function on a cut plane, 

a„2 log Z[X\ = - y dr^Gxir, f) . (92) 

The singular behavior of the variations of the Green function near the end points of the cut under the 
action of the symmetries gives us the mechanism to find an exact expression for logZ in terms of the 
solution of a second order non linear differential equation of the Painleve V type. This method was first 
used in pT] . 

Let us take the case A = e*^^" with a £ [0, 1). The Green function G{r,r') is uniquely defined by the 
following three requirements: 

a. - It satisfies the Helmholtz equation with a point like source 

{-Af + m^)G{r,r^)=S{r-r'). (93) 

b. - The boundary conditions are (they also hold for the Green function derivatives) 

lim G((a;,e),r ) = e'^""" lim G((x,-e)y ) for x e [L2, Li] , (94) 

Jim G(f,r') = 0. (95) 

|r| — *oo 

c- G{r,r') is bounded everywhere (including the cut) except at r = r*. 

This last requirement follows from the asymptotic condition (PO)) . We will write the Green function as 
G{z, z') as a shortcut of G(z, z, z' , z', Li, L2), where z and z are the complex coordinates x + iy and x — iy. 
It is Hermitian G{z,z') = G{z',z)*, and the time reflection symmetry gives G{z,z')* = G{z.z'). The 
reflexion operation 

R{x,y) = {Li + L2-x,y) (96) 
leaves the Helmholtz equation, the cut, and boundary conditions invariant. Thus we have 

G{z,z')^G{Rz,Rz'). (97) 

Due to the boundary conditions, near the end points of the interval [L2,Li] the Green function must 
have branch cut singularities. The requirement that the function must remain bounded on the cut and 
the equation imply that the most singular terms of G{z, z') for z near Li (and fixed z') have to be of 
the form 

G(z, z') ~ (z - Li)''^i(z') + (z - L,Y-''S2iz') . (98) 

Note that the contributions at this order must be analytic or anti-analytic in z in order to cancel the 
Laplacian term in (j93p . 

In the following the fact that, due to the uniqueness of the solution, a function which satisfies the 
Helmholtz equation and is bounded everywhere including the cut, must vanish identically, is used repeat- 
edly as a main argument in the calculation. This means that carefully analyzing the singular behavior at 
the extreme points of the cut of various quantities (formed by G(z,z'), iS'i(z), £'2(2), or their derivatives) 
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and combining them in order to cancel these singularities, one can actually construct equations which are 
valid everywhere. For example, from the analysis of the singular behavior of the derivatives Dl^G and 
SlsG at the singular points, using (|97p and ([55]) . we obtain the following fundamental relations 

(99) 
(100) 



dL,G{z,z') = S{z)^AS{z')., 
dL,S(z) = -fSiRz), 



where S is the vector with components Si and 5*2, 7 is a certain (unknow up to this point) matrix function 
oi L = Li — L2, and A is a constant Hermitian matrix. The solution for the half infinity cut which can 
be studied with standard methods gives A = — 47r(l — a)aai, with cri the Pauli matrix. It also holds for 
consistency of these equations, 



r 



Then, the equations 



and 



lead to the singular behavior 
S{z) 



1 

Itt 



(l-a)> 



(101) 



(102) 



for z in the vicinity of Li. 

The equations ([5^ and 
and 5*2 



with 



allows us to express the partition function in terms of the functions Si 
dLdm^ log Z = - y 5U 5 = Stto (1 - a) H{L) , (103) 



H{L) ^ / dxdy Sl{z)S2{z) . 



(104) 

Then, in order to compute logZ we need more information on S{z). 

With this aim, we exploit the symmetries that the Helmholtz equation has without imposing boundary 
conditions on the cut, to find relations for S and its derivatives. The idea is that due to translation 
and rotation symmetries, dyS and dgS, where dg — xdy — ydx is the rotation operator, also satisfy the 
Helmholtz equation and boundary conditions. Finding combinations which are free from divergences one 
finds the following equations 



dyS{z) 

dgS{z) 



i{7, asjSiRz) + ias,dxS{z) - £,S{z) , 
LidyS{z) + iqS(z) — iL^a^S{Rz) . 



Here ^ is another unknown matrix function of L and 



9 = 



1 - a 



(105) 
(106) 



(107) 



The R reflected equations to (|105p and (|106p also hold. The consistency of these equations with the 
Helmholtz equation and (|100p gives the algebraic equations 



£}A + A£^ ^ 

{{7, ^3}, a = 0, 



(108) 
(109) 
(110) 

(111) 



and the differential equation 



1 ^7 -1 
^cTs +7 97 + 9 + CT3 



(112) 



The algebraic equations (|10ip , (|108l - lllip for the matrices are solved using the parametrization 



TO 



u b 
c u 



if3i 
-ih 



(113) 
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Figure 5: The c„ functions for a scalar field, with, from top to bottom, n = 2, 3, and n — > oo. Here 
t = mL. 

where u, b, c, and (32 are real functions of t = mL, and + 1 = /3i/32- From this and the differential 
equation (jll2p it follows that all variables can be expressed in terms of u and u', and we have 

t 1 + M-^ + W"') 

This nonlinear second order ordinary differential equation can be transformed to take the form of a Painleve 
V equation [MllMj . 

In order to obtain a boundary condition for (|114p consider the Green function Go{z, z') of the Helmholtz 
equation without the cut, and the equation 

{Si{x)d^Go{x,xi) - Go{x,xi)d^Si{x)) - -5\x - xi)Si{x) . (115) 

Integrating this equation on the plane we have 

dx {Si{x)df,GQ{x,xi) - Go{x,xi)df_,Si{x)) = -Si{xi) , (116) 

where the integration contour is around the cut [L2,Li\. Then one can use a massless limit expansion 
Go(0, r) ^ — l/(27r)(log(rm/2) +7e), where is the Euler constant, and the massless limit solution for 
Si{x) of eqs. (|105m06|) . which is given in terms of hypergeometric functions, in eq. (|116p . in order to 
derive the boundary condition 

^ T^i T - - 1)* (logW + ""s) + ■■■ as t ^ , (117) 

t yiogt + Ks) 

with Ks = ^ log(2) + 27b + and il>[a\ the digamma function. u{t) admits a series expansion 

in terms of powers of t and \og{t) around the origin, with coefficients which are totally determined by Kg, 
and the differential equation. 

The long distance limit follows from the connection formulae for Painleve equations |55| . or the form 
factor expansion (see section 3.1.6) 

2 

Ua{t) ^ — ^yi^{oiT^)Ki^2a{t) ^S t — !■ OO . (118) 
TT 

So far, we have solved partially the problem since we still do not have an explicit expression for (|104p . 
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Figure 6: The function c{t) for a real scalar. At the origin it takes the value 1/3 = 0.33... and decays 
exponentially at infinity. The dashed curves are the leading short and long distance approximations. 



To find this quantity, we define the following auxiliary real integrals 



BiiL) = 


J dxdy Sl{z)Si{Rz) , 


(119) 


B2{L) = 


J dxdy S2{z)S2{Rz) , 


(120) 


Bi2{L) = 


J dxdy S2{z)Si{Rz) , 


(121) 


X,{L) = 


j dxdy Sl{z)Si{z) , 


(122) 


X2{L) = 


J dxdy 5*2(2;) 5*2 (2:) . 


(123) 



Using eqs. (|100l) . (|105p . and (|106p on these expressions, a complete set of liner equations among these 
capital letter variables can be found. As a result, all the capital letter variables can be written in terms of 
the functions u, u' by solving a linear system. In particular we get 

H={16TTa{l-a)my^tu^. (124) 

Combining p03p and (|124p we have 

/oo 
dyyu^y). (125) 

Therefore, the final expressions for the c functions (O and ([S]) are 

^ n— 1 



k=l 

00 



c{t) - -TT / db w^,,+,/2it), (127) 

2 Jq cosh(7ro)^ 

Eq. (|127p follows from ([^ by using A = e*^'^", and a ~ —ib+ 1/2, where b e (0,cxd). Wc have introduced 
a 1/2 factor in (jl26p and (jl27p because we want to present the results for a real scalar instead of a complex 
one. The function c„ is plotted in figure [5] for n = 2, 3 and n —^ 00. The function c{t) is shown in figure 
ini The leading long distance terms on the entropy functions are 

11 P p 

n — 1 47r t 
c{t) - ^tKi{2t). (129) 
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Note from p28p that the limits n ^ 1 and t — > oo do not commute (see further details in [301 [19]). Formula 
(|129p follows directly from (jllSp and (|127p . The short distance expansions read 

c„(t) ^ ^ + ^ + ^Oog^^W)^ (130) 

at) - (131) 

The constant term corresponds to the conformal case for which there is a general result identifying c = 
Cv/3, where Cy is the Virasoro central charge [551 [TU|. Here it is Cy — 1 since we are considering a real 
scalar field. 

The sharp cusp of c{t) at the origin, due to the l/log(i) term in (|13ip . is related to an infrared 
divergence for the entropy in the massless limit. In fact, for to this behavior gives for the entropy 
differences 

S{L) - 5(io) - dl ^ ^ - log(L/Lo) + - log (- log(m L)) - - log (- log(m Lq)) . (132) 

This suggests an infrared divergence S{L) ^ 1/2 log(— log(m)) for the entropy of any set. This fact can 
be checked by numerical simulations on the lattice, and has the following heuristic explanation. In the 
massless limit the homogeneous component of the field is a zero mode. The correlation function diverges 
logarithmically with the mass, and thus the typical size of the fluctuations on the homogeneous mode 
grows as (— log(TO))^/^. Correspondingly, the entropy grows as the logarithm of this volume in field space 
[56] . and becomes infrared divergent S{L) ^ l/21og(— log(TO)). This term in the entanglement entropy is 
independent of the number of components of the set, due to its infrared origin. In consequence, the mutual 
information I{A, S) ~ ^ log(— log(TO)) is also infrared divergent. 

3.1.2 Single interval for a massive Dirac field 

In this section we present a derivation of the entropic c function for a massive Dirac field by relating the 
Dirac Green function to the scalar one studied in 3.1.1. This derivation is different from the one in [30j . 
which is discussed later in 3.1.6. 

The functional Z[e^'^'^°'] can be calculated exploiting the relation between the free energy and the Green 
function 

dlogZ 

— ; =trG_D, (133) 

am 

where m is the field mass. We will take advantage of what we learned from the scalar case in the previous 
section. The Euclidean Green function Gd satisfies the equation 

(Tm^m + m)GD {z, z') = M(x - x')d{y -y'). (134) 

For definiteness we choose the Euclidean gamma matrices as 7^ = and 7^ = a^. For the scalar Green 
function Gs{z,z') we have 

(A - m^)Gs{z, z') = -5{x - x')8{y - y') . (135) 

Then 

CA ni'^)(rr.(7 7'] J m r.(7 7'))- {^''^'"'"l^'^y^y"! S'{x~x')S{y-y')\ , . 

(A - TO )(Gz5(z, z ) - 11 TO Gs{z, z )) - _ _ _ _ , (136) 

where we have used that 

(A - to^)Gd(z,z') = (5^7^ - TO)(a^7^ + to)Gd(z,z') 

= -f^S'{x - x')5{y - y') + -i^5{x - x')S'{y - y') - m6{x - x')5{y - y') . (137) 

We define the function T as the difference 

Gd{z, z')-Im Gs{z, z') = T{z, z') . (138) 

The sum of the diagonal components and the difference of the off diagonal ones in T satisfy the homogeneous 
Helmholtz equation according to (jl36p . Thus, they are non zero only if they are not bounded at Li or L2- 



19 




Figure 7: The Cn functions for a Dirac field and, from top to bottom, n — 2, 5, and n — > oo. 

According to (PT|1 the type of divergence which is possible for the fermionic Green function is of the same 
type than the ones of Si{z) and Si{Rz) defined in section 3.1.1. Then we can expand 

r(z, z') = aip{z, z') + a^qiz, z') + lw\z)f{z') + a2w\z)g{z') (139) 
where p, q, w, /, g are unknown functions, and w{z) ~ ( <^/e?\ ) ■ Having reached to this anzats, we can 



^SiiRz) 

now impose the hermiticity condition j^gI^j^ = Gd, where 7^ — ^7^7^, and the Dirac equation (jl34p for 
Gd = ^iTi Gs + T. This gives enough equations to determine the unknown functions in terms of quantities 
related to the scalar Green function. Making use of the equations found for the scalar case, we arrive at 
the following result for the trace 

trGc 2mtrG5 - 167ra(l- a)^ . (140) 

Pi 

Writing this in terms of u and u' and using (|133p we have 

with t — mL. Comparing this expression with the corresponding one for the scalar (jl25l) . it suggests that 
we could write it simply as 

/•oo 

Wa{t) = - / yv\y)dy, (142) 



where 

{l-2a)u{t)+tu'{t) 

" ^) = 77n — 27T\\ ' 143) 

t(l + u^(t)) 

and search for a differential equation for v. Surprisingly, this variable satisfies a similar differential equation 
of Painleve V type as the one satisfied by it, 

,, 1 , w ,9 -J 4a^ V 

v" + -v' = -- ^v'^ + -:r-. 9- 144 

t 1 — t'^ 1 — 

This follows from (|143p and the equation (|114p for u. The transformation l|143p . which maps functions 
satisfying Painleve equations to functions satisfying different Painleve equations is one particular case of 
the transformation group studied in [57j . The boundary conditions follows from the ones for m, 

2 

v{t) — > — sm{Tia)K2a{t) as t 00 , (145) 

TT 

v{t) -2a{\ogt+ kd) as t^Q, (146) 

where kd = -log(2) + 27E + ^^M+p^. 
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Figure 8: The c function for a Dirac field. Also shown with dashed lines are the long distance leading term 
and the series expansions around the origin up to orders t^, f^, t^, t^, and t^°. 



Thus, (|142p . (|144p and (|145p give the exact value of c„(t) and c{t), 

(n-l)/2 

1 - n 



c„(i) - -— J2 ^k/n{t). (147) 



fc=-(ri-l)/2 



oo 



Jq smh(7r6)^ 

The eq. (|148|) follows from using A = e*^'"' and a = ~ib, with 6 G (0, oo). The fmictions c„(t) are 
shown in figure [7] for some values of n. 

The leading long distance terms on these functions read 

Cn{t) -— +0( — ), (149) 

n — 1 ZTT t 

c{t) ^ \tKi{2t). (150) 

Eq. (|150p gives twice the corresponding term for a real scalar. Remarkably, it has been recently shown 
using the form factor approach, that the large distance leading term behavior c(L) ~ 1/4 miLKi{2'miL), 
where the sum is over the different particles of the theory counted with their multiplicity, also holds for 
the integrable interacting theories [33] , and also out of integrability [SS] . 

In order to improve the short distance expansion we can expand Va{t) by a direct use of the differential 
equations close to the conformal limit. We have the series solution of (I144[) around the origin 

C30 2m+l 

7;,(i)= ^ /™.„log"(0. (151) 

m— n— 

The full expansion requires only the knowledge of the constant term /o,o — kd- The first terms read 
Va{t) = ~2a\og{t) + KD +t'^ (^^{2a- 8a^ + KD - Sa'^KD - 4:aK% - K%) + (152) 

+ \ (-a + 8a^ + Sa^KD + Sukd) \og{t) + (-4a^ - So^kd) log^(i) + 2a^ log^(t)) + 0{t^) . 

The integration constant for Wa is given by the conformal limit Wa{0) = —2a}. The short distance expansion 
of the c- function follows from this expansion through p42p and (jl48p . It has the general form 

oo 2m 

c(t)= ^t^™^c™,„log"(0, (153) 
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Figure 9: Numerical vs. analytical results for the scalar (the lower curve) and Dirac (top curve) c-functions. 

where the Cm,n are functions of the The first terms are 

. i - + (1 + -'-«-"^;^'°«'^' ) - 0,163494 ,^ , (154) 

The expansion up to order t^" is plotted in figure [5] The series does not seem to converge for t > 1. For 
Cn the expansion starts with the leading terms 

Cn{t) ~ ^ - logf + Oif log(t)) . (155) 

3.1.3 Summary for one interval: Dirac vs. scalar fields 

Interestingly enough, when it comes to the entanglement entropy, the Dirac and the scalar c functions 
depend on one and the same differential equation. In order to see this, we only have to change variables. 
We summarize the results for a complex scalar and a Dirac field in the following way 

/•OO 











sinh(7r&)^ 
cosh(7r6)2 



= / db ^ ^^^ ws , (157) 



where 



W{D,S) = J dy y U(D^s){y) , (158) 

(D,S) {D,S) 

The differences between the fermionic and scalar cases are exclusively due to the boundary conditions 



2 

Unit) — smh{bTT)Ki2b{t) as t ^ oo (160) 

TT 

Unit) -2b{\ogt+ kd) as t^O (161) 

2 

us{t) ^ — cosh{bTT)Ki2b{t) as t ^ OO (162) 

TT 

usit) ^ ~\ asi-.0, (163) 

t(logi + Ks) 

where kd ^ - log(2) + 2^e + and ks ^ - log(2) + 2^e + '/'[V2+i^+V-[i/2-.^ ^ 

We show in figure[5]a comparison of the exact analytical c-functions for a scalar and Dirac fields, along 
with the results from the lattice numerical simulations (see section 3.1.8). The bosonic c-function has 
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Figure 10: The plane cut along the intervals {ui, Vi) with i — 1, ..p 



the same central charge at the conformal point as the Dirac field, and thus they both tend to 1/3 at the 
origin. The scalar function then rapidly goes to half the Dirac one (that is, approaches the c functions 
corresponding to a Majorana field). Note also the very different behavior at the origin. The bosonic 
c- function has a 1/ log(i) term that can be ascribed to the zero mode which is present at the conformal 
point. 

3.1.4 Entropy saturation at long distance 

At large t — Lm the function c{t) = LdS{L)/dL quickly goes to zero, and the entropy stops to grow. For 
this saturation limit one has 

f°° c(v) 

S(oo) = / dy — = -c(0) log(me) + const . (164) 

Jem y 

This result about the dependence on log(m) of the saturation constant for the entropy was established with 
all generality in [TU]. For the scalar and the Dirac field we have c(0) = 1/3. For the scalar, a subleading 
term ^ log(— log(m)) term is also present. 

3.1.5 Dirac field: bosonization and the massless multicomponent case 

Consider a massless Dirac field and a general set V consisting in a collection of disjoint intervals {ui,Vi), 
i = l,---,p (figure [TU]). In order to calculate S{V) and Sn{V) we follow section 2.1.2. Accordingly, the 
problem is reduced to n decoupled fields living on a single plane. These fields are multivalued, since 
when encircling Cm or C^. they are multiplied by e'"^'^ or e~*"^'^, respectively. 

That multivaluedness can be disposed of, at the expense of coupling singled-valued fields 'if'' to an 
external gauge field which is a pure gauge everywhere, except at the points Ui and Vi where it is vortex- 
like. Thus we arrive to the Lagrangian density 

£k = (^M + « K) + • (165) 

The reverse step would be to get rid of the gauge field by performing a singular gauge transformation 

(where xq is an arbitrary fixed point). Since the transformation is singular, one goes back to a multivalued 
field. In order to reproduce the boundary conditions on ^I''^, we should have 

/ dx^A^ix)^-—, i dx^Alix)^ — . (167) 

Equations in (jl67|) hold for any two circuits C^. and Cy. containing Ui and Vi respectively. Thus 

k ^ 

e^^d^Alix) = 27r- ^ [5{x - u{) - 5{x - v^)] , (168) 
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where the presence of a vortex- antivortex pair for each k and each interval is expUcit. Then, the partition 
function can be obtained as vacuum expectation values in the free Dirac theory 



where is the Dirac current, A'^ satisfies (|168L and we adopted a normalization such that (1) — 1. 

In order to evaluate (|169p . it is quite convenient to use the bosonization technique in two dimensions 
[30j . One expresses the fermionic current in terms of a dual scalar field as 

jj^^J-e^'^d,^. (170) 

Then, the functional becomes 



where the vacuum expectation values correspond to the scalar theory. For the massless fermion this is 
simply a free massless scalar 

1 

2' 



C4, = -d^^d^(j). (172) 



Since is quadratic 



with the correlator 



G{x^y) = -^\og\x-y\. (174) 
It follows that (|17ip can be written as 

logZ[e*2.fe/„] ^ ^S(y), (175) 

^(^) = ^log|Mj - -^log|u, - Mjl -^log|wj - -ploge. (176) 

i,i i<j i<j 

Here e is a cutoff introduced to split the coincidence points, \ui — Ui\, \vi — Vi\ e. Summing over k and 
using ([H]) and ([7]) we obtain 

Sn = ^^S(l/), (177) 

S = ^EiV). (178) 

This formula first appeared in [lOj as a proposal for the entanglement entropy of disjoint intervals for any 
conformal theory in two dimensions. It was later discovered that this is not the case However, it 

is the correct formula for the free massless fermion model. It is interesting that (|178p does not seem to 
hold for the a model dual to the free fermion given in terms of a compactified scalar 59 . The reason for 
this is probably that while the theories can be mapped to each other, this mapping is non-local, and the 
definition of what is the algebra of operators in region differ between these models (60) . 

Equation (|178p has an interesting corollary: recalling the definition ^ for the mutual information, it 
follows that, for non-intersecting sets A, B and C 

I{A,B^C) = I{A,B)+I{A,C). (179) 

That is, in contrast to the entropy, the mutual information is extensive for the massless free fermion. This 
curious property does not hold in the massive case or in more dimensions |61) . 
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3.1.6 Scalar and Dirac partition functions as Sine-Gordon correlators 

We can also use the bosonization technique to deal with (|169p for the partition function for the case of a 
massive fermion. This gives place to an alternative way of obtaining the partition functions already studied 
in sections 3.1.1 and 3.1.2. As in the massless case, we still have eq. (jl7ip for the partition function. Now, 
however, the bosonization of the massive fermion theory leads to a sine-Gordon theory at the free fermion 
point [55], with Lagrangian 

= \ {pi.4>d''4> + A cos( \/4^0)) , (180) 

where A is a dimensionful parameter. 

This allows us to identify the partition functions with sinc-Gordon correlators of vertex operators. In 
the case of a single interval we have 

Z[e*2.T ~ {VaV-a) , (181) 

where Va is : exp[j\/47ra(/)] : and a G (—1/2,1/2). These type of correlators have been studied in the 
literature, starting from the one for a = 1/2 which is related to the Ising model spin and magnetization 
correlators The correlator (jl8ip has a long distance expansion in terms of intermediate particle states, 
called form factor expansion. This reads 

/. g.V4¥a^W .. g.V4¥aV(0) A = du,...du2n 

n=Q ('^■) 

x/a(Ml, ■•■,U2n)/a'(w2n, ■■■,"!) , (182) 

where fa{ui, ...,U2n) is the form factor 

/a(^l,...,.2„) = (-l)"("-l)/2 (^^)" (n i^y) ^ (183) 

and 

A(Mi,...,U2n) = 2„ - ' . (184) 

Y{r=lY\s=n+l {Ur + Us) 

The first term on this expansion gives place to the long distance behavior pSOp for the entropy [301 133j ■ 
It can be shown that this series may be summed as a Fredhlom determinant and satisfies the Painleve 
eqs. (|14mi46|) [M]. This proof follows a very different route than the one in sections 3.1.1 and 3.1.2. The 
form factor expansion and integrable field theory techniques have been successfully used to approximate 
the entanglement entropy in one dimensional integrable models [33l [35l [19] . 

For a scalar field we also have an expression of the partition function in terms of sine-Gordon correlators, 

^ I , (185) 

with a g (0, 1). This can be checked directly using the Painleve equations, and psip . It also follows from 
the form factor expansion for the scalar case [65]. 

3.1.7 Dirac field: small mass expansion in the real time approach 

In order to compute the entropy for the multicomponent massive case, we find convenient to use the real 
time approach. First, we express the real time approach entropy formula (j75p in terms of the resolvent 
i?= (C- 1/2-1-/3)-! as 



S{V) =- dptr 

Jl/2 



(/3-l/2)(i?(/3)-i?(-/3)) 



/3+1/2 

The correlator for a two dimensional Dirac field reads 



(186) 



1 1T1 1 m 

C{x, y) = -5{x ~y) + —K„{m\x ~ y\)^^ + — K,{m{x ~ y)h' , (187) 

Z ZVr ZTT 
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Where 7" and 7^ are the Dirac matrices, and 7^ = 7*^7^. In the massless hmit this gives 

Co{x,y)^^S{x^y)l + ^^^'. (188) 

Fortunately, the resolvent for this integral operator inside a region formed by n disjoint intervals (ui,Vi) 
is known from the theory of singular integral equations . We have 

= _ 1/4) Uf3S{x-y)+^- -— j, (189) 

where 

z(x) = logf-nki^^) . (190) 

With this resolvent at hand we are in position to compute the entropy for the massless case and to make 
expansions for the massive case |51j . 

In (jl86p the term proportional to the identity cancels with the corresponding one in the resolvent (jl89p . 
Then we have, 



Sm=o{y) — / dl3 dx lim 



1/2 



sm 



^log(f+^) iz{x)~z{y)) 



V y-^- iP + 1/2) (x - y) 



(191) 



Integrating over /3 first, this gives 

iM_£Mcoth((z(a;)-z(y))/2)-l 1 



f^^ mi^Z eoth((.(.) - ziy))/2) 1 I f y ( 1 
Jv y^"" - v){z{x) - z{y)) Gjy 



Uj X — Vi 



\ ( X! h - - ^ log \u, -uj\-Y^ log \v, -Vj\- nloge 



(192) 



This corresponds to the eq. p77|) obtained with the Euclidean method [50] , 

The short distance expansion of the correlator up to second order in the field mass is C = Co + Ci + 
C2 + ... where 

C.{x,y) = -£(7.+log(!lfc^))7°, (193) 
C.{x,y) = '-^{x^y)U-UloJ^!^\],^. (194) 



47r ' \ 2 

The perturbative expansion of the resolvent for small mass is 

Rv{l3) = Rl{l3) - R"y{f3)CiR'^{l3) - i?°.(/3)C2i??,(/3) + i?^(/3)Cii?^(/3)Cii??,(/3) - ... . (195) 
Using this in p86|) we get (see [FT] for details) 

S{V) = S„.^oiV) - ^ log^(m) + log(m) P^^^^'^ " " ' m^L^, (196) 



6 o V y . OV . ^ g 

-2m^ J dx J dy log |a; - y\ 6(z{x) + z{y))^ + 0{m^ \og{mf) , 

where Lt = ~ ""O the total length of the intervals. For a single interval this coincides with the 

corresponding terms in the expansion (|154p obtained with the help of the differential equation. The last 
term is the first one which gives a non-extensive mutual information j61j . 
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3.1.8 Lattice Hamiltonian and correlators in one dimension 

In this section we describe the Hamiltonian and correlators which can be used in a one dimensional lattice 
to calculate numerically the entanglement entropy. The numerical results of the c functions for a scalar 
and a Dirac field are compared with the analytical ones in figure O showing a perfect accord. 
Scalar 

We take the lattice Hamiltonian for a real massive scalar as 

n— — OO 

We have set the lattice spacing to one. The correlators (I55l) and are 

7^ „ix{rn—n) 

dx , , (198) 

TT 47r^m2 + 2(1 - cos(a:)) 



f da;— e"('"-"V"^^ + 2(1 



cos(x)) . (199) 



Dirac fermions 

We take the Hamiltonian for a Dirac fcrmion discretized on the lattice 



oo / ■ \ 

^ = E {^h"l\^n+i - *„) - h.c.) + mv&h"*™] • (200) 



n— — oo 



The two dimensional matrices 7" and 7^ can be taken at will while satisfying the fundamental relations 
for the Dirac matrices (i.e. 7'^ = cti and 7^ ~ i(72)- The correlator (|67p is 



2 



47r-y/ + sin^(a;) 



Due to the fermion doubling on the lattice, one has to divide the lattice results by 2 in two dimensions in 
order to get the entropy corresponding to a Dirac field in the continuum limit. 



3.2 Two dimensions: logarithmic term for polygonal sets 

In the continuum limit described by a quantum field theory the entanglement entropy is divergent due to 
the presence of an unbounded number of local degrees of freedom. The singularities structure is encoded in 
the expansion ([3]) where the dimensionless coefficient g^iV) of the logarithmic term is particularly relevant 
since it is the only one universal. Logarithmic divergent terms in the entropy have been previously found 
in four dimensional black hole space-times [23^ . They are present generically in even spacetime dimensions 
{d odd) for sets with smooth curved boundaries. This follows from the heat kernel expansion for conical 
manifolds with smooth singularity surface (see section 3.3.3). In [6T and [5S| it was shown that there 
is also a logarithmic term in d = 2 for sets V with non-smooth boundary (see also [HI [25]). The figure 
[TT] shows the logarithmic term for a square in a two dimensional lattice. Since go{V) is dimensionless, 
extensive and local on the boundary, for V a spatial polygonal set, it must be of the form 

ffo(^)=E«(^^)' (202) 

where the sum is over all vertices Vi and Xi is the vertex angle. On general grounds one also expects 
point-like vertex induced logarithmic terms in any dimensions. 
The alpha-entropies also contain a logarithmic term 

^"liog = E^"(^*)l°g(^^)' (203) 

analogous to the one in the entropy, with s{x) — lim„_^i s„(a;). In eq. (|203p A is a parameter with the 
dimensions of an energy, depending on V and on the particular theory. For a massless field it is the inverse 
of any typical dimension i? of ^ and when the mass dominates, MR ^ 1, it can be taken as A = M . 

In this section we review the analytic results on s{x) obtained in [S7] and [55] for a free scalar and 
Dirac fields respectively, and show the results match with the numerical simulations on a two dimensional 
lattice. 
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Figure 11: The logarithmic term for the entropy of a massless scalar, where F is a square of side L. Here 
we have plotted the results of the numerical evaluations of the entanglement entropy on the lattice, where 
we have subtracted the linear and constant terms in a fit of the form 5* = cq + ciL — slogL. The solid 
curve is —0.0472 logL, showing perfect accord with the numerical data. 



3.2.1 Scalar field 

In order to calculate the logarithmic coefficient s„ for a free scalar we consider a plane angular sector as the 
simplest set with vertex contributions. We start studying the Green function as in 3.1.1. The eigenfunctions 
of the Laplacian admit separation of variables, and by direct calculation of the radial component of the 
eigenfunctions, we are left with a two dimensional reduced problem. This later consists in finding the trace 
of the Green function on a sphere with a cut with particular boundary conditions on it inherited from 
the original problem. The Myers method, used to find the entropy of a scalar field in a plane with a cut 
in 3.1.1 is well suited to the present case as well. The study of the divergences structure of the Green 
function in the vicinity of the set boundary and the identification of the symmetries are again the basic 
ingredients. 

Dimensional reduction 

We start considering the Green function G{fi,r2) for a complex scalar of mass M in three Euclidean 
dimensions subject to the boundary conditions ITTl) . To be explicit, we have (see figure [T^ 

+M2)G(ri,r2) = J(ri - rl) , (204) 
lim G(ri +e7],r2) = e*^^" lim G(ri - ef), fa) , ri^V, (205) 

e^0+ £—0+ 

where fj is orthogonal to the plane of V . This is related to the functional Z according to 

5M2logZ[e''2™] = - y" dr3G(f, r) . (206) 

The Laplacian and the boundary conditions allow the separation of angular and radial equations in 
polar coordinates. Using standard methods we arrive at the expression 

x-^ /■ A Ji+y(Ari) Ji , ^(Ara) 

G{n,r,) = Y, J dX^^^M0i,Vim02,^2) ^^"^^ J-fJ^ \ (207) 

where J is the Bessel function. Here the sum is over the normalized eigenvectors iJji,{9,Lp) of the angular 
equation 

An^'u ^ + , (208) 

where Ao is the Laplacian on the sphere, Ao = ji^^(sin0^) + sin^ e ' ^i^h domain given by the 
functions satisfying the boundary conditions inherited from (|205l) . The precise expressions for and v 
are not relevant in what follows. 
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Figure 12: Sphere with a cut at the intersection with the plane angular sector of angle x. 



Taking the trace j dr^G{'r,r) in eq. (|207p gives (disregarding an unimportant divergent constant) 



dM^ log^ = - ^ + 1/2) = ^tr^-Ao + I . (209) 

V 

Though this expression is divergent, the piece we are interested in, which is the one dependent on the 
angle a;, is finite. To proceed, it is convenient to express the trace of the square root of the operator in 
(j209p in terms of the corresponding resolvent. We have the identity 



1 1 ( 1 
tr^/-Ao + - = - / A^tr- ^ -d\. (210) 



Green function on a sphere with a cut 

The problem is now reduced to the calculation of the trace of the two dimensional Green function on a 
sphere with a cut of angle a;, where the boundary conditions inherited from (j205p are imposed. 

This problem is the analogous on the sphere to the one solved for the plane in 3.1.2. Following the 
Myers method step by step as in 3.1.2, we find the analytic expression for the trace of the Green function 
as a solution of a system of ordinary differential equations. The details of the derivation are given in [67j . 
Explicitly we find 

trG^^) = 87r(l - a)a \ H{y)dy . (211) 



Here the function H{x) is the solution of the following set of ordinary non linear differential equations 
(we omit the subscript a and the dependence on x of the variables for notational convenience) 

771 

H' = -- {bB2 + cBi+2uBi2) , (212) 
X[ = -m (6Bi2 + uBi) , (213) 

X'2 = -771 {CB12+UB2) , (214) 

c' = — 2 m /32 u csc(a;) sin(x/2) — c (1 — a) csc(a::) (1 + cos(a::)) , (215) 
h' = — 2 m/3i u csc(a;) sin(a::/2) — 6a csc(a::) (1 + cos(a::)) , (216) 

771 1 
u' = -y sec(x/2) (6/32 + c/3i) + -u tan (j;/2) , (217) 

where Bi, B2, B12, f3i, P2 are functions of x given in terms of Xi, X2, c, 6, and u by the following set 
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Figure 13: Coefficient of the vertex induced logaritlimic term in the entropy for a scalar field as a function 
of the vertex angle x. The three points shown correspond to numerical evaluations on the lattice for 
X = 7r/4, -it/2, and 37r/4. 



of algebraic equations 

cos(x/2) 
87ra(l — a) 
sm{x/2) 



?7ra(l — a) 



sin(.T/2) H -m{(3iX2+ P2X1) + 2m cos(a;/2) u B12 , 



- cos(j:/2) H ~m tan(a;/2) (/Ji^a + /32^i) + m sin(x/2)(6B2 + cBi) , 



(218) 



(219) 



= -msin(a;/2)(cXi - 6X2) + m tan(a;/2)(/32Si - f3iB2) + (1 - 2a) cos(x/2)Bi2 ,(220) 
= -4a(a- l)-m2(4-8/3i/32 + 6c + 3u2) 

—4 cos(x) (a(a- 1) + 771^(^2 + 1)) +m^cos(2x)(6c- u^) , (221) 
= (2a- l)ucos(a;/2) +mtan(a;/2)(/3ic- 6/32) ■ (222) 



The boundary conditions at a; — s- tt are 



H{7: 

b{Tr 
c{tt 



= 0, 



r(-a) (cosh (^) Im (1 + a + f )] - f sinh (^)) 
22"/i (cos (2a7r) + cosh(7r/^)) r(l + a) |r (i - a + f^) f 

-'^1(7'') I a^(l-a) , 



= 0, 



)l-2o. 



,(l-a)|r(i + a+f)| 



mV'^{\ + a) 

K'^) I a-»(l-a) , 



(223) 
(224) 

(225) 
(226) 

(227) 
(228) 



where ^ = V 4 nrfi — 1 and t/" is the digamma function. The meaning of the extra variables -Bi, i?2, S12, 
Xi, X2, u, 6, c, /9i and /32 is the same as in 3.1.2. The trace in (j21ip is regularized such that it vanishes 
when cc = TT, where there is no vertex point and no logarithmic term is present in the entropies. 

Gathering aU the results together, using eqs. (^05]) with A = M, HH), ((^Dg)) . (PIU)) . and (I^TTI) . we 
arrive at the final result for the logarithmic coefficients 



8k{n-k) 



and 



s(a;) 



^1 n2 (n - 1) Ji/; 
167r(t2 + 1/4) 



rfm rn (m^ — 1/4) 2 / dy Hk_{y,m) . 



cosh (tt^) Ji/2 



dm TO (j 



1/4) 



(229) 



(230) 



where the function _ff is solution of the above set of ordinary non linear differential equations. 
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In [70j the partition function for a Dirac fermion on the Poincare disk with multiphcative boundary 
conditions imposed on a geodesic segment has been written in terms of a solution of the Painleve VI 
differential equation. Though we were not able to find an explicit relation, it is possible that our results 
for H{x) could have also an expression in terms of solutions of these type of equations. The function H 
also gives the exact entropy functions for a massive scalar and a spatial segment in 1 + 1 dimensional de 
Sitter space [67] . 

An economic way to integrate the equations is to expand the functions involved in (|212j|228"| in Taylor 
series around x — tt, and obtain analytically the coefficients using the differential equations. Then the 
above integrals over t and the mass can be done for each coefficient separately. With this method we 
have produced the curve of figure [T51 which show ss{x) up to order {x — tt)^^. Some of the coefficients 
are tabulated in table [TJ In the picture are also plotted the values of ss for x — 7r/4, tt/2 and 3/Att 
obtained by numerical simulations in the lattice. They show a perfect accord (around one percent error) 
with the analytical results. These particular values of the angle are the ones for which the coefficient 
can be calculated with very small error on a square lattice of limited size (in the present case it was 
200 X 200 points). The numerical methods consist of evaluating the entropy for a massless Dirac field (see 
[57J) for a given shape (square, triangle, etc.) and different overall size A, and then fitting the result as 
S ^Co + CiX + C-i A-i + C_2 - ss log(A). 

The small angle limit of s{x) relates the problem on the cut sphere with the corresponding one in the 
plane, treated in 3.1.1. One has from this mapping 

, , f„°° dt c{t) , , 

six)^^-^ 231 

where c is the entropic c-function for a scalar. This result can be understood in more general terms. This 
is explained in section 3.3.1. 

3.2.2 Dirac field 

The problem for a Dirac field and Y a plane angular sector can also be dimensionally reduced by separation 
of variables. Then, as in the scalar case, we have to calculate the trace of the Green function on a two 
dimensional sphere with a cut. This problem is then related to the scalar one already treated in 3.2.1. 
The partition function for a Dirac field ^ in three dimensions is 

^[e^2.aj ^ y p^tp^g- / dr' *ti,3* ^ (232) 
where I?3 is the Dirac operator given by 

= {I'd, + ii) , (233) 

and 7* = cr' are the Pauli matrices and /i the mass. The boundary condition for the spinors is 

*+(r) = e''2'^"*-(r) , r (^V . (234) 

Here 5*+ and are the limit values of the field on each of the sides the two dimensional angular sector 
Y has in three dimensions. 

The functional Z is calculated exploiting the relation between the free energy and the Green function, 

^ = TrGg) , I?3Gg)(r, ^) = S'{r- r^) . (235) 

Dimensional reduction 

The Dirac operator and the boundary conditions allow the separation of the angular and radial equations 
in polar coordinates. We use this fact in order to reduce the problem to one in two dimensions. In this 
coordinates writes 

V3 = -+ rdr + M , D= [fde + td4>) , (236) 

r 

and where the redefined gamma matrices 7 are 

^^=r^a\ ^^^r^a\, = l^l^^ = sin07V • (237) 

OXi OXi OXi 
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Figure 14: Coefficient of the vertex induced logarithmic term in the entropy for a Dirac field as a function 
of the vertex angle x. The three points shown correspond to numerical evaluations on the lattice for 
X = 7r/4, -it/2, and Sn/A. 



The dimensional reduction proceeds as in the scalar case, by integrating out the radial variables in the 
Green function. For details see [68] • We arrive at 

trCg) =-i-tr|r^-l/2| . (238) 

Using (|235p and (|236p the logarithmic divergent contribution to the partition function is 

logZliog = -\tr\YD - l/2|log(e/i). (239) 

The trace in (|239p can be calculated using the integral representation in terms of the resolvent [BS] 

tr|7''D- 1/2| = / dmmir{i{YD-l/2) + m)-K (240) 

The operator 

X>2 = i{YD - 1/2) + m (241) 

is a two dimensional Dirac operator on the sphere where the parameter m plays the role of a mass. 

Thus, we have to find the trace of the Green function of a Dirac field on a two-dimensional sphere with 
a cut. This satisfies 

V.G'^^^ ^ ^ 5\z - z') , (242) 

and the boundary conditions (|234p . 

Relation between the scalar and Dirac Green functions 

As in the flat space analog of sections 3.1.1 and 3.1.2, we can relate the Dirac Green function on the sphere 
to the Green function Gg of a scalar field of mass M. Defining an auxiliary quantity G as 

G{z,z')^vlG^p{z,z'), (243) 

we have that 

V2Giz,z') = ^ 5\z-z'), (244) 
where the scalar and the fermion masses are related by 

= i + . (245) 

Thus, the difference 

G^^\z,z')^G{z,z') = Q{z,z') (246) 
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''2 






J-) 




s(V2) 


s(3/47r) 


ss 


7.8125310-^ 


5.45402 lO-"* 


5.3465610"'' 


5.40167 10"'^ 


7.9410-^ 


0.02366 


0.005040 


sd 


7.8125310-3 


5.0142610"'* 


4.8129910"^ 


4.8552310-6 


7.2210-2 


0.02329 


0.005022 



Table 1: The first four non zero Taylor coefficients of so{x) and ss{x) (a complex scalar) for x around tt, 

^(^) ~ S'^j'^'* ~ ^'^'^ the coefficient of the term 1/x of these functions for x ^ 0, s{x) ^ c'^l/x. 
The value of the functions for x — tt/2 and x = 3/4 tt are also shown. The quadratic coefficients of s/)(x) 
and ss{x) coincide. 



satisfies the Dirac, 'DQ{z,z') = 0, and Helmholtz equation, (— + AP)Q{z, z') — 0, without sources. 
Therefore it would be identically zero if it where bounded. Q{z, z') is however unbounded at the extreme 
points Li and L2 of the cut. These singularities can be conveniently eliminated by linear combination with 
the functions Si{z). Therefore expanding the most general for of Q{z,z'), imposing the Dirac equation, 
and using the differential equations for Si{z) and S2{z) given in [67j (eqs. (81-84) of that work), one 
derives the Dirac Green function in terms of quantities related to the scalar one. 

Here we need only the part of the trace of G)j = G + Q which is odd in the mass to, since the terms 



even in to do not contribute to the integral (|240p . We have 

2m tr G^s' ^ — \ , " ,' ' ^, ^ . (247) 



tr {x, TO, a) 



(2) 167ra(l - a)TO (4/3iXi cos{x/2) - b Bisin^ix)) 
odd "•^ M {Af3i^ - 62sin2(a;)) 



Here tr G5 , b, /3i, Xi and Bi arc functions of x which are the analogous on the sphere to the quantities 
defined on section 3.1.1 for the plane. They are given in terms of the system of algebraic and ordinary 
differential equations (j218ll222"|) . 

The resuh follows combining ([^55)1 and We have 

sd{x) = I dt ^ , ^ / dmm trG^^\x,m,-it) . (248) 



2 sinh-^(7ri) 

The relevant part of tr Gd{x, to, —it) is given by (|247p and the formulae at eqs. (|212ll228p . where we have 
to make the replacement a — —it. The imaginary part cancels in (|247p as it should. 

The function sd is plotted in figure [HI with some lattice results. Table [T] shows the coefficients of the 
Taylor expansion around x — t: for ss and sd. Also shown is the coefficient of the term ^ 1/x for s{x) in 
the limit a; ~ 0. In this limit the formula (j23ip applies, with c{t) the c-function corresponding to a Dirac 
field in one dimension. 

3.2.3 Lattice Hamiltonian and correlators in two dimensions 
Scalar 

In the lattice simulations we have taken the following lattice Hamiltonian for a real scalar in three dimen- 
sions 

n.m— — CO 

We have set the lattice spacing to one. The correlators (1551) and ([5^ are 



(249) 



^ n J J cos(ia;) cos(jy) 



/ dy ^ ' . , (250) 



87r2 ^2(1 - cos(a;)) + 2(1 - cos(y) + m^) 



1 



(tto.otTij) — ^— J / dx j dycos{ix) cos{jy)\/2{l — cos{x)) + 2(1 — cos(y)) -|- to^ . (251) 
Dirac fermion 

We take the lattice Hamiltonian for a Dirac fermion as 
^ = ((*L,„7%'(^™+i,„ - *m,„) + *L,„7V(*™,„+i - *„.,„)) - h.C.)+m¥^^„j'>^m,n ■ (252) 

71,771 
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The correlators are 



(*u.^..> = + r d. r dy !!ni±pM2^±^mi,^i^i^-,)^yi'^-^^^ . (253) 

J-^ J-^ 87r2 Jto2 +sin^(a;) +sin^(y) 

Due to the fermion doubhng on the lattice, one has to divide the lattice results by 4 in two dimensions in 
order to get the entropy corresponding to a Dirac field in the continuum limit. 

It is relevant to the accuracy of the entropy calculation to evaluate the correlators with enough precision. 
This can be very time-consuming. We have found it is much faster to evaluate one of the two integrals in 
the correlators analytically (in terms of polynomials times elliptic functions) using a program for analytic 
mathematical manipulations, and then doing the last integral numerically. 



3.3 Universal terms in the entanglement entropy for higher dimensions 
3.3.1 Dimensional reduction 

For free fields, some universal terms in the entanglement entropy in higher dimensions can be obtained from 
the low dimensional results ^5SJ. Let us consider sets in fc -I- d spatial dimensions of the form Z — B x X , 
where i? is a box on the first k coordinates xi, ...Xk, of sides having lengths Ri, i — I, fc, and X is a set 
in d dimensions. We are interested in the entropy of Z in the limit of large Ri, which should be extensive in 
the sides Ri. Thus, we can compactify the directions Xi, with i — 1, fc, by imposing periodic boundary 
conditions Xi = Xi + Ri, without changing the result of the leading extensive term. For a free field of mass 
m we Fourier decompose it into the corresponding field modes in the compact directions. The problem 
then reduces to a d dimensional one with an infinity tower of massive fields. For example, for a free scalar 
we obtain the tower of fields 4'qi,---,qk^ where qi,...,qk are integers, such that the momentum in the first fc 
components is pi — ^Qi- From the point of view of the non compact Xk+i, Xk+d directions, these fields 
have masses given by 



k 



27r 



=1 

Summing over the contributions of all the decoupled d dimensional fields we have 



S{Z)^ ^ S{X,M,,_,,). (255) 

{gi} = -oo 



In the limit of large Ri we can convert these sums into integrals 

kA 



^ 2Wf(i72TT) / dpp^^'SiX^V^^^')^ (256) 







where the transversal area A — JliLi ^i- The universal terms in S{Z) will then come from the ones of X 
after this integration over the mass. For non scalar fields the factors of the spin multiplicities have to be 



incorporated into (j256|) . 

One important application of this formula is the case of a thin set Z which is large in all directions 
excepting one, with size L. The leading universal term in the entropy is proportional to the area A in the 
large face transversal to L, 

6S{Z) = -K^^. (257) 

The same coefficient k gives the dominant part of the mutual information of two sets A and B with large 
parallel faces of area A, at a distance L much smaller than the other dimensions, 

I{A,B)^K^. (258) 

In order to calculate k we use (|256p where X is a one dimensional set of length L. In consequence 
LdS{Z)/dL is related to the entropic c-function of the one dimensional case. It follows that 



2d-2^^r(^) 



t Jo 



<t) ^ . ^-l^.._l^ / dyi yr" / dy2 yt" c Jyl + , (259) 
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d 


«;(bosons) 


K(fermions) 


2 


3.9710-^ 


3.6110-^ 


3 


5.5410-3 


5.3810-3 


4 


1.3110-3 


1.3010-3 


5 


4.08 10-"^ 


4.0610-'' 



Table 2: The coefficients k for massless bosons and fermions (one field degree of freedom) for different 
spatial dimension d. 



where t = mL, and m is the mass of the field in d spatial dimensions. Here the function c{t) = (riBCBit) + 
npcp(t)), where ub and np are the multiplicity of the bosonic and ferniionic degrees of freedom, and 
CB{t) and cpit) are the one dimensional entropic functions corresponding to a real scalar (boson) and a 
Majorana fermion fields (see section 3.1.3). In the massless limit k reduces to a constant 

n^({d~l) 2'-' n'^rf^)) ' Tdy c{y) . (260) 







The coefficient k can be calculated numerically with good precision. Table (j3.3.1|) shows k for the 
first few dimensions. For large dimensions the bosonic and ferniionic coefficients approach each other, 
what is a consequence of the fact the c-function has the same asymptotic large t behavior in both cases, 
c{t) ~ jtKi{2t). Using this approximation in (I260p we have for large d 

— ^ 2 ; (261) 



This has less than the 3% error already for d = 5. 

For d ~ 2, and considering now V as an angular sector of small angle x instead of a long strip, one 
should add the contribution of the term (|257p along the length of the angular sector, 

6S ^ -K [ dL^ ^ -loge. (262) 
J xL X 

For d = 2 eq. (|260p gives k = , and (|262p gives the logarithmic term in the entropy (|23ip for the 

small angle limit ^67]. 

3.3.2 Two component sets in the large distance limit 

Another interesting universal term is the correction for the entropy of two sets A and B which are far 
apart |61j . This is measured by the mutual information I{A,B), which goes to zero for long separating 
distance L. We find convenient to use the real time approach in order to compute this limit. For a free 
Dirac field we can use the formula (jl86|) for the entropy in terms of the resolvent i?(/3) = (C + /? — ^) 
of the field correlator C. For a massless field we have 

C(x,y) 0|vl'(a;)vl/(y)t|0 h{x ~ y) + ^c7V /^_~fj^, , (263) 

z \x y\ 

where c is a constant which depends on d, and are the Dirac matrices. 

Considering V ^ AU B with a large separation vector L between A and B, the resolvent can be 
expanded perturbatively, 

RviP) = R°v{f3) - R'^(j3)C,R°y{p) + R'^{l3)C,R"y{p)CiR''y{(3) - ... . (264) 

Here Ry{j3) is the unperturbed resolvent, 

RUP)={''i^^ V (265) 



Rb{P) 

and Ci is the field correlator evaluated for the separation vector L 



y/L- f Ia.b ^ 
=«c^7^ . . . . (266) 
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Figure 15: The log-log plot of the mutual information I{A, B) for two squares A and B of side 10 lat- 
tice points as a function of the separation distance L in two dimensions. The line with larger slope 
(J(^, B) ^ L^*) corresponds to a Dirac field, while the other to a scalar one {I{A, B) ~ 



The kernel 1a,b{x^ y) ~ 1 for any x ^ A and y £ B. 

The second term in (I195P does not contribute since it has zero trace. The third term is proportional 
to the squared of the field correlator and leads to 

I [A, B)^^ r dp {(3 - 1/2) [ tr (/3) + R^{f3)R^i(3)) - (/3 ^ -/?)] , (267) 



1/2 

where V means the set V after an inversion of coordinates followed by a reflection in the plane perpendicular 
to L. The bar over the resolvent and the square of the resolvent means sum over the spatial variables, 
Ox = ^^dx dyO{x,y). Thus, we have the power I{A,B) ^ L^^"^. The coefficient should have length 
dimension 2d, and depend on A and B and their relative position but not on the distance L. Moreover, it 
must be monotonically increasing with A and B. For the massive case I{A, B) is also proportional to the 
field correlator squared and decays exponentially quickly. In particular, in one dimension we have |51| 

IiA,B)^^m^LfLfiK^imL) + Kl{mL)) {1 + 0{La,b/L) + ...) , (268) 

where L\ and are the total lengths of A and B, assumed much smaller than and L. 

A similar calculation for scalars on a lattice shows also a decay going like the field correlator squared. 
For massless fields this gives I{A, B) ~ L^^^'^^^\ which is slower than the fermionic case. It is instructive 
to obtain this behavior with a different approach. The mutual information is an upper bound on correlation 

m 

/(A,i.)>iii^^^^l^i^^l^, (269) 

where Oa and Ob are any bounded operators in the local algebras corresponding to A and B, with norms 
||0^|| and II Ob II respectively. We cannot apply this relation directly to the smeared scalar field operator, 
/ dx a{x)4>{x), with a{x) an smooth function of compact support, since this is not bounded. However, using 
the inequality (j269|) with the unitary operators Oa = exp {i J dx a{x)(j){x)) and Ob = exp(i J dy f3{y)<p{x)), 
and a{x), P{x) vanishing outside A, B respectively, we obtain, 

liA, B) > cons ^0(O)<?!>(L)J>^ - L-^'^'^-^l (270) 
3.3.3 Logarithmic term for smooth boundaries in four dimensions 

In the Euclidean approach the entanglement entropy is obtained from the effective action W — — log(Z) 
calculated in an Euclidean d -1- 1 dimensional manifold with conical singularities on dV. We have the 
general expansion 

W = non universal + a^+i log e + finite . (271) 
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The non universal terms diverge as inverse powers of the cutoff e. fn the case of a free theory a^+i coincides 
with the corresponding heat kernel coefficient ([55]) . For a general conformal theory it is determined by 
the integrated conformal anomaly. In order to see this, consider an infinitesimal rescaling of the metric 
gf"" ^ (1 - 25\)g'"'. Since 

in terms of the energy momentum tensor T^i/, we have 

^ = -/'i^'V?(n)' (273) 

which is the integrated trace anomaly (classically, for a conformally invariant theory we should have 
T'^ = 0). On the other hand, due to the conformal invariance of the action, scaling the metric as above 
must give the same result as keeping the metric constant but contracting the cutoff as e (1 — 5X)e. Using 
(l?7T1) this gives 

ad+i = I dx' (T^^) . (274) 



The trace anomaly for smooth manifolds has a local expression in terms of polynomials of the curvature 
tensor of dimensionality equal to d + 1. Thus, it must vanish for odd spacetime dimensions. In four 
spacetime dimensions it writes '71' 

04 = —^{aE-icI), (275) 
57D07r^ 

E - J dx^R^f.^.R'^f^'''' - AR^.R^-' + R^) , (276) 

/ = Jdx^ (i?a;3M.i?"^^'' - 2i?^,i?^'' + ^ (277) 

where a and c are the anomaly coefficients, which are a characteristic property of the theory. Here they 
are normalized such that a = c = 1 for the real scalar. For a Weyl fermion it is a = 11/2, c = 3 j7l j . 

As they stand these equations cannot be used for evaluating the entanglement entropy in flat space 
since the bulk curvature tensor Rap^m = 0, and, on the other hand, there is a conical singularity at 
dV . However, it was argued by Fursaev and Solodukhin in 23] that at lowest order in the deficit angle 
5 = 27r(l — n), the contribution of the surface term is proportional to this deficit, 

a4~(n-l)s, (278) 

and s can be calculated using the same expressions (|275II277|) . with a curvature tensor proportional to a 
delta-function concentrated in dV . From ([7]) this gives the logarithmic term in the entanglement entropy 



'loe 



sloge. (279) 



The value of s in (|278p was obtained in [55] for the case of zero extrinsic curvature of the surface dV , 
and used to calculate corrections to the entanglement entropy of the black hole [531 [3S] (the extrinsic 
curvature of the horizon vanishes). The contribution of the extrinsic curvature is however crucial for the 
geometric entropy in flat space. This was obtained recently by Solodukhin in [40]. His approach consists in 
demanding conformal invariance of the most general geometric expression for the logarithmic contribution, 
and then using the holographic anzats of Ryu and Takayanagi (which is applicable to the special case of 
the entanglement entropy of SU{N) superconformal gauge theories [75]) in order to calibrate a coefficient. 
The result in flat space i^ [ID] 

{k^K, - + ^ / (K^K, - Jfcf^fc;,) ■ (280) 

dV Z4U7r Jgy z 



Here fc^^ = ~l'jllvda'n}p is the second fundamental form, with z = 1,2 are a pair of unit vectors 
orthogonal to dV ^ and 7^j/ — 5^^^ — ri^nj, is the induced metric on the surface. In the first term, the 



^The relation of the present normalization of the anomaly coefficients with the one in go] is A = a/(907r2), B = c/{30Tr^). 
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integral Jgyik'i^'^ kl^ — fcf^fc^^) = 47rx(9y), where x{9V) is the Euler number of the surface, which is 
purely topological. 

For a sphere ((MI)) gives [ID] 

^ = 1)' 

while for a cylinder of length L and radius R it is [JU] 

s = — -. (282) 
240 i? ^ 

The sphere and the cylinder are then sensitive to different anomaly coefficients. 

Note here that for the electromagnetic field it is a = 62 and c = 12. This shows the entanglement 
entropy of the electromagnetic field cannot be assimilated to the one of a pair of massless scalars. Also 
note that there are some conflicting results in the literature which do not coincide with the above formulas 
for the logarithmic coefficient in three spatial dimensions for scalar fields [73]. 



4 Outlook 

The obvious missing subject on this review about entanglement entropy for free fields is the one of free 
gauge fields. In fact, though there are a few works on this sense in the literature, the entanglement entropy 
for gauge fields has not been fully understood yet. There are some problems in the definition of S(V) in 
this case, since there is no gauge invariant partition of the Hilbert space as a tensor product of the Hilbert 
spaces corresponding to degree of freedom on two complementary regions (see for example [Zl])- Perhaps 
this is the origin of the boundary term noticed in the early works for the electromagnetic field [321 , which 
would make the area term non-gauge invariant. At present it is not known if the mutual information is 
well defined in the presence of gauge symmetries. 

Other problems directly related to free fields also seem to be worth of further work. One is the problem 
of defining the real time formulation for bosons directly in the continuum. In this case it is easy to see 
by dimensionality arguments that the correlators X and P of section 2.2.1 are not well defined in the 
continuum, while XP, which is the operator relevant to the entropy calculation, has the right units. One 
should be able to obtain this continuum operator directly, without passing through a regularization of X 
and P. We think this issue has to be understood before attacking the problem posed by the localization 
of Abelian gauge fields. 

Other interesting research subjects include the treatment of excitations, which seem to allow a system- 
atic formulation in the real time approach, and the development of a perturbative theory for the reduced 
density matrix. 
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